Secondary instability of electromagnetic ion-temperature-gradient modes for zonal flow generation 



Secondary instability of electromagnetic ion-temperature-gradient modes for zonal 
flow generation 

Johan Anderson, 1 ]^] Hans Nordman, 2 Rameswar Singh, 3 and Raghvendra Singh 3 

^Department of Nuclear Engineering, Chalmers University of Technology, 
SE-412 96 Goteborg, Sweden 

^Department of Earth and Space Sciences, Chalmers University of Technology, 
SE-412 96 Goteborg, Sweden 

3 ) Institute for Plasma Research, Bhat, Gandhinagar, Gujarat, 
India 382428 

(Dated: 19 January 2013) 

An analytical model for zonal flow generation by toroidal ion-temperature-gradient 
(1TG) modes, including finite ft electromagnetic effects, is derived. The derivation is 
based on a fluid model for ions and electrons and takes into account both linear and 
nonlinear /3 effects. The influence of finite plasma on the zonal flow growth rate 
{Izf) scaling is investigated for typical tokamak plasma parameters. The results show 
the importance of the zonal flows close to marginal stability where ^zf/iitg >> 1 
is obtained. In this region the parameter ^zf/iitg increases with j3, indicating that 
the ITG turbulence and associated transport would decrease with (3 at a faster rate 
than expected from a purely linear or quasi-linear analysis. 

PACS numbers: 52.30.-q, 52.35.Ra, 52.55.Fa 

Keywords: Zonal Flows, Ion- Temperature-Gradient Modes, Finite /3, Transport 



) Electronic mail: anderson.johan@gmail.com. 



1 



Secondary instability of electromagnetic ion-temperature-gradient modes for zonal flow generation 
I. INTRODUCTION 

The dependence of energy transport on plasma (3 ((3 = 8irn (T e + Ti)/B 2 , kinetic-to- 
magnetic pressure ratio) is of major importance for the operation and performance of a 
magnetic confinement fusion device. For tokamaks, the combination of high confinement 
and high would allow for a high fraction of bootstrap current as well as high fusion gain 
and offer the promise of a more compact and economically feasible tokamak reactor operating 
in steady-state. The investigation of such advanced confinement regimes is presently a high 
priority research area in both experimental and theoretical fusion plasma physics. 

Experimentally, the scaling of confinement with (3 has shown inconsistent results. In the 
commonly used empirical scaling law IPB98(y,2), a strong degradation of confinement with 
increasing (3 is predicted 1 . In dimensionless scaling experiments on the other hand, where 
(3 is varied while the other dimensionless variables are kept fixed, the strong degradation of 
confinement has not always been confirmed^. The mixed results may be due to the different 
turbulent types/regimes in the edge and core plasmas. 

Theoretically, it is well known that the interplay between ion temperature gradient (ITG) 
mode turbulence and zonal flows play a crucial role for the level of turbulent transport in 
core plasmas. Zonal flows are radially localized flows that are driven by the turbulence 
and propagating mainly in the poloidal direction. The zonal flows provide a strong shear 
stabilization of the turbulent eddies and are hence important for the self-regulation and 
saturation of the turbulence. 

Theoretical studies of ITG turbulence and transport including finite (3 electromagnetic 
effects have usually been based on linear and quasi-linear (QL) theories. It is known that 
the linear ITG mode growth rate is reduced by electromagnetic fluctuations, resulting in 
a favourable scaling of confinement with f3 in QL theories^— . Less studied is the role of 
electromagnetic effects for the generation of zonal flows. Nonlinear simulations of ITG turbu- 
lence including electromagnetic effects and zonal flow dynamics have been performed using 
both gyrokinetioii~— and gyro-fluid models^ 1 ^. Recently, nonlinear gyrokinetic simulations 
of ITG turbulence has reported a significant reduction of transport levels with increasing j3 
which could not be explained by the linear physics aloneii^. 
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From a theoretical point of view there are several analytical models for treating multi- 
scale interactions. Among the widely used models are the coherent mode coupling method 
(CMC), the wave kinetic equation (WKE) approach and the reductive perturbation expan- 
sion method. In comparison, the CMC model is based on a finite number of test waves, 
such as pump waves, zonal flows and side bands whereas the WKE analysis is based on 
the coupling of the micro-scale turbulence with the zonal flow through the WKE under the 
assumptions that there is a large separation of scales in space and time.- 



In the present paper, electromagnetic ft effects on I 



ysed based on a two-fluid model for the ions (Refs. ll9M2ll) and the electrons employing the 



G turbulence and transport is anal- 



WKE model for the zonal flow generation. Note that the kinetic ballooning mode (KBM) 
is however not included. A system of equations is derived which describes the coupling 
between the background ITG turbulence, using a wave-kinetic equation, and the zonal flow 
mode driven by Reynolds and Maxwell stress forces. The work extends a previous study 



(Ref. |22| ) by self-consistently including linear as well as nonlinear ft effects in the derivation. 

The derived dispersion relation for the zonal flow is solved numerically and the scalings of 
the zonal flow growth rate with plasma parameters, in particular with plasma ft, is studied 
and the implications for ITG driven transport scaling with ft are discussed. 

The paper is organised as follows. In Sec. II the fluid model used to describe the 
electromagnetic ITG modes are presented. The derivation of the zonal flow growth rate in 
the presence of a background of ITG modes is described in Sec. III. In Sec. IV the results 
are presented and discussed. Finally the conclusions are given in Sec. V. 



II. ELECTROMAGNETIC TOROIDAL 
ION-TEMPERATURE-GRADIENT DRIVEN MODES 



We will start by presenting the ion part of the fluid description used for toroidal ion- 
temperature-gradient (ITG) driven modes consisting of the ion continuity and ion temper- 



ature equations by following the Refs 



10 and 



22 



The electromagnetic effects enter through 



electron fluid model via quasi-neutrality while the ion branch is identical to the electrostatic 
case. Combining the ion and electron fluid model through quasi-neutrality results in the 
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dispersion relation for ITG modified by finite (3 effects. It has been found that the effect of 
parallel ion motion is weak on the reactive ITG modes and therefore it is neglected, moreover 
the effects of electron trapping is neglected for simplicity. The linearized ion-temperature 
and ion-continuity equations can be written 

uh = (u Di - u+i)4> + uj Di (n + fi) - k]_(u - 0^)0, (I) 

~ 2 ~ 2 _ ~ 5 ~ 

ioTi = {-u Di - riiU±i)<fr + -UDi(n + T t ) + -u Di Ti - k]_(u - 0:^)0. (2) 

Here we have assumed r = Tj/T e = 1 and ft = (L n / p s )5n/n , <p = [L n j p s )e8(f>/T e , % = 
(L n / p s )5T i /T i0 are the normalized ion particle density, the electrostatic potential and the 
ion temperature, respectively. We have denoted the ion diamagnetic and magnetic drift 
frequency as = —u+ e and upi = —oJoe where the geometrical quantities are calculated 
using a semi-local model where [ff) is defined by ujd {9) = u*€ n gi (9), with cj* = k y v± = 
kyp s c s / L n , see equations f TT7|) - ffTO]) below. The normalized gradient scale lengths are defined 
as Lf = — (dlnf /dr)^ 1 (/ = {n, Tj}), rji = L ra /L^, e n = 2L n /R where R is the major radius 
and CKj = (1 + T]i). The perpendicular length scale and time are normalized by p s and L n /c s , 
respectively. Here p s = c s /Vt ci where c s = ^jT e /rrii and Vt ci = eB/rriiC. We will start by 
deriving the linear ion density response of the form hi = Q(p for this system of equations. 
Combining Equations (fTT) - (}2"]) and eliminating the temperature perturbations we find a 
relation between the ion density and potential perturbations, 
T ~ 

iii = —(j) = Q<p, (3) 
10 5 

N = u 2 - —u Di uj + -u%t, (4) 

f 7 5 5 \ 

T = iu(u Di - + u+iUj Di (- - T)i) - -u 2 Di - (u + aiU*i)(uj - -u Di )k]_ \ . (5) 

Now we turn our attention to the electron fluid model. We will consider a low-/? tokamak 
equilibrium with Shafranov shifted circular magnetic surfaces while omitting the parallel 
magnetic perturbations (the compressional Alfven mode) and we will make use of a electric 
field representation of the form, 

£ =- v ^-cir e '" < 6 > 

where </> is the scalar potential, Am is the parallel component of the vector potential and eii is 
the unit vector along B. We find from the parallel momentum equation for electrons while 
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neglecting electron inertia a relation between the electron density, potential and parallel 
vector potential, 



ck\\ 

We will now use the quasi-neutrality condition (hi = h e ) in combination with the parallel 
electron momentum ([7]) and the ion density response ([3]) to determine the parallel vector 
potential in terms of the electrostatic potential yielding 

in =ch\ 1 - - $ = tt a (p. (8) 

00 - 0J± P 



Here, 



n a = cku- — (9) 

CO - 



In order to arrive at the dispersion relation we need 
electrostatic potential and the vector potential, see Ref. 
equation to find 



yet another equation relating the 



10 



We use the electron continuity 



um e = (u* e - u De )4> + u De P e - —k\\J\\, (10) 

ne 

where uo e is the electron magnetic drift frequency and P e = nT e + T e h e is the linearized 
electron pressure perturbation. Furthermore we assume that the electron parallel heat con- 
ductivity is large ViiT e m where Vy is taken along the total magnetic field line giving 

fc||T e = 77e^i||. (11) 

In the regime vu « ^ ~ Va « v te the parallel current is primarily carried by electrons 
(va = ^t^ nm . an d v t are the Alven and thermal speed, respectively) resulting in, 

fell J|| = fcy J|| e = ne\ [u* e - u))k\\(p H A\\ I (12) 

The parallel current density (Jm) can be determined through the parallel component of 
Ampere's law, 

Air 

VlA {l = --J h (13) 
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yielding the second relation between the potentials and A\\ as 

k \\4 = 1 + 1 7 S ~ A W 14 

Combining Equations (jHJ) and ffblj) and normalizing with a; 2 e we find the ITG mode dispersion 
relation as 

- k\k\ + (5 (±\ f(Q -1)(Q- e n9i ) + Ve e n g t - = 0. (15) 

Here, we have normalized the ITG mode real frequency and growth rate as u = The 
geometrical quantities will be determined using a semi-local analysis by assuming an ap- 
proximate eigenfunction while averageing the geometry dependent quantities along the field 
line. The form of the eigenfunction is assumed to be 10 , 

W(0) = + cos#) with |0|<7r. (16) 

V 37T 

In the dispersion relation we will replace fcn = (fcii), = (k±) and wo = (a;^)) by the 
averages defined through the integrals, 

(kl) = £ dBWj! = k 2 e (l + ^ (vr 2 - 7.5) - ^sa + ^a 2 ) , (17) 

/ n /2 5 5 \ 

d9^u D y = e n co, I - + -s - —a \ , (18) 

(k {l kl *„) = £ de*fc|,*ifc B * = (l + s 2 (y - 0.5) - 8 -sa + . (19) 

Here a = (3q 2 R{l + r/ e + (1 + r/j)) /(2L n ) and /3 = 87rn (T e + Ti)/B 2 is the plasma /3, g is 
the safety factor and s = rq'/q is the magnetic shear. The a-dependent term above (in Eq. 
[T71) represents the effects of Shafranov shift. We will now study the non-linear generation of 
zonal flow induced by toroidal ITG modes modified by electromagnetic effects. 



III. THE MODEL FOR ZONAL FLOW GENERATION 

In order to determine the zonal flow generation from the non-linear coupling of ITG 
modes modified by electromagnetic effects we will need to describe the temporal evolution 
of the zonal flow through the vorticity equation, 

(V- J) = (V- Jl) + (V- J\\) (20) 
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Here (f(k x ,k y )) = J d 2 kf(k x , k y ). The vorticity equation consists of two parts, the first 
including a derivative perpendicular to the field line and the second along the field line. At 
first, we will consider the two contributions separately where the perpendicular part can be 
written, 

(V ■ J_l) = (V ■ {eriiVu - rieVte)) + (V • en^) 
d\7 2 , $ 

= -^-([0,V ± 0]> (21) 
Here \A, B] = ir^^P — ^-W- is the Poisson bracket. Next, we consider the contribution from 

l ' J ox ay ay ox 

the parallel derivative of the current density, 

,1 , 1 1 en x Viii 

— VJy = — V|| Jj| L -fT JL ■ VJ II 

en en en B 

'^) 2 <[4vi-4,]> 

^ (id-o) (-^-) i 2 [0,vi« 

= (\Q a \ 2 [4>,Vl4>]) (22) 

where we have used equation (jHJ) to substitute the vector potential by the electric potential 
and we have asssumed that the variation along the field line is very small V||oJ|| ~ 0. The 
time evolution of the zonal flow potential ($) is them given by, 

|vi$ = -((l-|fi Q | 2 )[0,Vi0]) (23) 

We will now compute an estimate for the generation of zonal flows through the Reynolds 
stress and Maxwell stress terms. We consider the Reynolds stress, 

ft Vi0]> = -V^<|^> = -V 2 x Re j d 2 kk x k y \4>\ 2 . (24) 

Here Re stands for the real part of the integral and the gradient in X acts on the spatial 
scale of the zonal flow. We will now assume that there exist a wave action invariant of the 
form \<f)\ 2 = C(k x , k y )Nk- Now the zonal flow evolution becomes, 

^Vi$ = -V 2 x Re j d 2 kk x k y (l - \n a \ 2 )C(k x , k y )N k . (25) 

In order to close the system of equations we need an additional relation for the action 
invariant (Nk) which is given by the wave kinetic equation. The wave kinetic equation (see 
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Refs U, |22|, 1244291) for the generalized wave action N k in the presence of mean plasma flow 
due to the interaction between mean flow and small scale fluctuations is 

d d ( t x dN k (x,t) d (% _ \ dN k {x,t) 

= 7i TG N k (x, t) - AuN k (x, tf (26) 

In this analysis it is assumed that the RHS is approximately zero (stationary turbulence). 
The role of non-linear interactions among the ITG fluctuations (here represented by a 
non-linear frequency shift Aw) is to balance the linear growth rate. In the case when 
JlTGNk(x,t) — AojN k (x,t) 2 = 0, the expansion of equation ([26"]) is made under the assump- 
tion of small deviations from the equilibrium spectrum function; N k = N% + N k where N k 
evolves at the zonal flow time and space scale (Q, q x , q y = 0) of the form ^ ~ e n^x-m for 
^ = {$,N k }, as 

d - dN° 
-i(n q -q x v gx + i 1ITG )N k = g^{k-zx V(l + Vi)$)-^. (27) 

While using Equation ( 1271) we now find the perturbed action density as, 

Nk ~ (n q - qx v gx r + 1 * ITG k « ( v * (1 + W (28) 

and substituting equation (128|) into the zonal flow evolution we obtain, 

^V 2 $ - -V 2 Re fd 2 kk k ^ l - W)Vl{l + Vl) ^ C(k k (29) 
g V ± $- V x Rejdkk x k y {Qq _ qxVgx) 2 + 1 2 TG C{k x ,k y ){ Qk Y (29) 

The adiabatic invariant iV° = j s determined by the energy density E k and the real fre- 
quency i>j rk . An approximate wave action density can be obtained using the same methodol- 
where the linear electron equation (TH1) and quasi-neutrality flE]) are used to 



ogy as in Ref . 



22 



find the modified normal coordinates if = + a k Ti for which a generalized invariant is found 
of the form N k = |<£>| 2 = C(k)\(p\ 2 . It is generally found that C(k) is only weakly dependent 
on the wavevector k as long as the FLR effects are small. The remaining integral displayed 
in Equation (1251 can be solved in the two limits < 1 and ^^^^ > 1. We assume a certain 

i x v 7/TG 1ITG 

spectral form on the action density N® and that C(k x ,k y ) is weakly dependent on k x . In 
the limit where the linear growth rate is much larger that the zonal flow growth < 1 
we find the dispersion relation, 

7/tg ok x 



, . _ k x kl{\ - \Vt a \ 2 ) dN? 
V q = u£(l ~ ql) d 2 k »\ ' aU C(k x , k y ){-°-±). (30) 
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Here the 1 — \fl a \ 2 represents the electromagnetic effects on the zonal flow evolution and the 
f — q 2 term is the FLR stabilization. We choose the particular form of the the spectrum as 
a Gaussian wave packet in k x with width A and delta function in k y such that, 

N° k ^N e-^ s)1 6(k y -k y0 ). (31) 

We have chosen a drift wave packet centered around the most unstable mode in k y and a 
spectrum in k x similar to that used in Ref. [30|. Now the derivative on the action density 



is easily found and the integral can be computed resulting in the final zonal flow dispersion 
relation as, 

n q = i{K - \n a \ 2 )ql{l - ^)S^|0| 2 . (32) 

Ik 

Here the dispersion relation may be modified with the inclusion of the FLR non-linearities 



(see e.g. Ref. 



221 ) by setting K = 1 + r + t5 p where 



_ A,A 2 + | 7 | rG 



Ai = (r]i ~ I) + I^itg, (34) 
5 

A 2 = UJlTG + ^ e n9 (35) 

and ujitg is the real frequency and 7/tg is the linear growth rate of the ITG mode. 

Next we will consider the more interesting limit > 1 in the integral (129|) . In this 
limit the zonal flows are expected to have an impact on the ITG turbulence. We assume 
that the coefficient C(k x , k y ) is weakly dependent on k and that the group velocity can be 
written v gx = k x k y f(rji, r) e , (3, . . .) where f(rji, r) e ,(3, . . .) is independent of of k) similar to the 



case in Ref 



22 



We can now rewrite the integral as 

*n,u u ^ / , 2 , >^(1-I^| 2 ) ( dNl 
(fig - q x v gx ) dk x ' 



n q = q 2 x (l-q 2 x )C(k Xl k y ) j d 2 k- 



,2,1 ,.2, C(k x ,k y ) f 2 k y v gx (l \Q a \ 2 ) dN®, 



f(r]i, r] ei p,...)J {il q - q x v gx ) dk x 

We consider the same spectral form as in Eq. EH and performing one partial integration the 
dispersion relation is readily found, 

(n g - q x v gx ) 2 = -q 2 x (K - \n a \ 2 )(l - q*)k y AV^\4>\ 2 (37) 
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In the following we will use a fixed width of the spectrum with Aa/7t = 1.0 corresponding to 
a monochromatic wave packet in k x . In the following section we will explore this dispersion 
relation numerically and discuss the results and its implications. 



IV. RESULTS AND DISCUSSION 

The algebraic equation ( |37|) describing the zonal flow growth rate is solved numerically 
with the ITG eigenvalues ujtg taken from a numerical solution of the ITG dispersion relation 
( [15]) . The zonal flow growth rate is then normalized to the ITG growth rate to highlight the 
competition between the linear ITG drive and the stabilizing effect of the zonal flow mode 
through shearing of the turbulent eddies. The results are expected to give an indication of 
the strength of the shearing rate u s ~ (where <£> is the zonal zonal flow component of 
the electrostatic potential) relative to the linear growth rate. The results for 'Jzf/iitg is 
calculated for a turbulence saturation level, corresponding to a mixing length estimate with 
|r = y~t~ = ^° ^ na ^ * s fi xe d-— In experimental tokamak plasmas the core density profiles are 
rather flat (e„ > 1) whereas the edge profiles are peaked (e n < 1) with a typical experimental 
value of the plasma around 1%. 

In Fig. 1 the ITG eigenvalues (normalized to the electron diamagnetic drift frequency) 
function of are displayed with parameter for e n = 2.0, T] e = 0, r = 1, q = 1.5, 

s = 0.5 and k\p 2 = 0.1. The results are shown for rji = 3.5 (dashed line), and for rji = 4.0 
(solid line). In the electrostatic limit — > 0, the analytical results as obtained from Eqs 8 



- 10 of Ref. 



22 



(neglecting FLR effects) is u/uj+e = —3.03 + i0.90 for the case rji = 3.5 and 



uj/u^e = — 3.03 -Ml. 28 for rji = 4.0, in good agreement with the numerical results of Fig. la. 



The results show that the 
previous studies (J8 



12 



and 



TG growth rate is reduced with increasing as expected from 



21|). 

The corresponding results for the zonal flow growth rate (normalized to the ITG growth 
rate) versus f3 are shown in Fig. lb. As observed, the normalized zonal flow growth increases 
for large close to marginal stability. The scaling illustrates the competition between the 
linear and nonlinear stabilizing effects. Close to marginal stability (for rji = 3.5), the 
decrease of 7/tg due to increasing dominates resulting in a normalized zonal flow growth 
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rate that increases with (3. This would indicate that the ITG turbulence and transport 
decrease faster with (3 than expected from a purely linear analysis, in agreement with recent 



simulation results (11 ll). For larger rji, a decreasing zonal flow drive is observed due to Maxwell 



stress. We note that this result differs with that reported in Ref. (|31| and |32|). using the 
drift-Alfven wave branch neglecting effects of curvature, where the zonal flow growth initially 
decreases with /3 and reaches a minimum and then increases. The consequence of such a /3 
dependence is not apparent, but could indicate a transition to a more favourable confinement 
regime for (3 > f3 min . 



Figure 2 illustrates the ITG and zonal flow growth rates as a function of rji with f3 as 
a parameter. The other parameters are e n — 2, rj e — 0, r = 1, q — 1.5, s = 0.5 and 
kj_p 2 = 0.1. The results are shown for (3 = 0.1% (dash-dotted line), [3 = 0.5% (dotted 
line), (3 = 1.0% (dashed line) and (3 = 1.3% (solid line). The linear ITG growth rates 
(normalized to the electron diamagnetic drift frequency) in Fig. 2a illustrates the typical j3 
stabilization with a /9-threshold close to the analytical result rjuh = 3.01 for (3-^-0 (from 
Eq. 10 of Ref. Q). Fig. 2b displays the corresponding zonal flow growth rates (normalized 
to the ITG growth rate). The un-normalized zonal flow growth is weakly dependent on rji, 
resulting in a normalized zonal flow growth rate which strongly increases when approaching 
the stability threshold rjuh- The results indicate the importance of the zonal flows close to 
marginal stability (rji < 4.0) where Jzf/iitg >> 1 is obtained. In this region the effects 
of zonal flow increases with increasing (3 whereas for larger rji the opposite trend is found. 
This is in line with the strong nonlinear upshift of the critical ion temperature gradient 
with increasing (3 and converging Dimits shift for larger (3 recently observed in nonlinear 
gyrokinetic simulations of ITG turbulence^ 1 ^. However a complete treatment of the Dimits 
shift requires a model for the zonal flow saturation mechanisms which is outside the scope 
of the present paper. We note that the same trend, with an increase of ^zf/iitg with 
increasing (3 is obtained in Figure (lb) close to marginal stability (r]i = 3.5). For larger 
r/i > 4, the condition ^zf/iitg > 1 is however not satisfied. 

11 
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V. CONCLUSION 

A system of fluid equations describing the coupling between the zonal flow mode and the 
background ITG turbulence including finite electromagnetic effects is derived. The model 
equations include the linear stabilization of the ITG mode due to finite electromagnetic 
perturbations as well as the nonlinear effects on the zonal flow entering through the 
Maxwell stress force. The scaling of the zonal flow growth rate with plasma parameters is 
studied and the implications for ITG driven transport are discussed. It is found that the 
ZF growth rate relative to the ITG growth rate increases with close to marginal ITG 
mode stability. The result indicates a stabilization of the ITG turbulence and transport 
at a faster rate than expected from a purely linear or quasi-linear analysis. Such behaviour 
has recently been observed in nonlinear gyrokinetic simulations of ITG turbulence^ 1 ^. The 
results are also line with the strong nonlinear upshift of the critical ion temperature gradient 
with increasing and converging Dimits shift for larger recently observed in nonlinear 
gyrokinetic simulations of ITG turbulence^. We note that close to marginal stability the 
increase of ^zf/iitg is dominated by the linear stabilization of 7/tg whereas for larger 
rji a decreasing zonal flow drive is observed due to the competition between Reynolds and 
Maxwell stresses. 

In the immediate future it is of interest to complement the present model by including 
geometry effects and A± which is important in high plasmas relevant for spherical systems. 
In addition, further investigations of zonal flow stability and saturation mechanisms and their 
relation to transport barriers are left for future study. 
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FIG. la. Numerical solution to Eq. (I15jl displaying the ITG growth rate (normalized with w^ e ) 
versus f3 for r\ e = 0, q = 1.5, s = 0.5, k\p 2 = 0.1, Ay^ = 1.0 with e n as a parameter. Results are 
shown for rji = 3.5 (dashed line) and rji = 4.0 (solid line). 
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FIG. lb. Zonal flow growth rate (normalized with 7jtg) versus j3 with r\i as a parameter for the 
same parameters as in Fig. la as obtained numerically by solving Eq. (|37|) . A fixed ITG turbulence 
saturation level <f> = <f>o was used. 
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FIG. 2a. Numerical solution to Eq. (|15p showing the ITG growth rate (normalized with u;* e ) versus 
rji for rj e = 0, e n = 2, q = 1.5, s = 0.5, k\p 2 = 0.1, = 1.0 with /3 as a parameter. Results 

are shown for j3 = 0.1% (dash-dotted line), (3 = 0.5% (dashed line), j3 = 1.0% (dashed line) and 
= 1.3% (solid line) 
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FIG. 2b. Zonal flow growth rate (normalized with 7jtg) versus r\i with j3 as a parameter for the 
same parameters as in Fig. 2a as obtained numerically by solving Eq. (|37|) . A ITG turbulence 
saturation level (ft = <f>o was used. 
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